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The  ability  to  predict  the  aerodynamic  characteristics  of  a  flight 
vehicle  without  actually  fabricating  and  testing  the  vehicle  has  always 
been  a  highly  sought-after  engineering  goal.  Original  predictions 
generally  involved  numerous,  simplifying  assumptions  and  were  done  via  hand 
calculations.  Today,  however,  the  computer  has  greatly  expanded  the 
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with  the  modern  capabilities  of  a  high-speed  computing  machine  in  order  to 
investigate  a  topic  of  very  current  interest.  This  topic  is  the  ejector 
wing:  and  the  goal  has  been  to  develop  a  methodology  that  can  be  used  to 
calculate  the  two-dimensional  lift  per  unit  span  of  any  ejector  wing 
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Abstract 


The  theoretical  framework  and  general  solution  procedure  have  been 
developed  for  calculation  of  the  lift  per  unit  span  of  an  ejector  wing 
model.  This  model  is  based  on  the  fundamentals  of  Potential  Flow  and  Thin 
Airfoil  Theory  and  incorporates  a  point  sink,  point  source,  and  two  bound 
vortex  sheets.  The  solution  procedure  consistently  satisfies  both  the  flow 
tangency  and  Kutta  conditions,  but  its  usefulness  is  shown  to  be  highly 
dependent  on  the  number  of  control  points  used.  Numerical  examples  are 
presented  for  cases  involving  five  control  points,  and  the  lift 
calculations  which  result  are  shown  to  be  inconsistent.  A  FORTRAN  computer 
program  is  included  for  the  five  control  point  case,  but  it  can  be  modified 
to  accomnodate  any  number  of  control  points. 
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THEORETICAL  DETERMINATION  OF  THE 
LIFT  OF  A  SIMULATED 
EJECTOR  WING 


I.  Introduction 


Background 

An  ejector  wing  looks  much  like  a  conventional  subsonic  aircraft  wing 
except  for  the  presence  of  one  or  more  slots  which  run  along  its  span-wise 
direction.  Figure  1  shows  an  actual  ejector  wing  that  has  recently  been 

/ 

tested  in  the  7'  x  10'  Army  tunnel  at  the  NASA  Ames  Research  Center  (7). 

The  figure  shows  the  underside  and  leading  edge  of  the  wing  and  the 
presence  of  four  distinct  slots.  These  slots  occupy  what  would  be  the 
interior  region  of  a  conventional  wing,  and  they  exit  at  the  top  surface 
near  the  trailing  edge.  Figure  2  is  a  cut-away  drawing  of  the  wing  showing 
both  its  external  and  internal  configurations  (7).  Figure  3  shows 
schematic  representations  of  the  cross  sections  of  both  conventional  and 
ejector  wings. 

The  principle  upon  which  the  design  of  an  ejector  wing  is  based 
involves  the  injection  of  high  energy  air  into  the  slots.  This  air  is 

I 

obtained  from  the  aircraft's  exhaust,  and  it  is  mixed  with  the  air  flowing 
through  the  slots.  When  this  mixture  reaches  the  wing's  upper  surface,  it 
acts  to  delay  the  onset  of  boundary  layer  separation.  By  keeping  the  flow 
attached  to  the  ejector  wing's  surface  over  a  greater  area,  an  increase  in 
lift  over  that  of  a  conventional  wing  is  obtained. 

The  particular  wing  design  shown  in  Figures  1  and  2  was  first  proposed 
by  the  Flight  Dynamics  Laboratory  of  the  Air  Force  Wright  Aeronautical 
Laboratories  at  Wri ght-Patterson  Air  Force  Base,  Ohio.  It  was  designed 
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Figure  1.  Leading  Edge  and  Ift^erside  of  an  Ejector  Wing 


Ejec 


«- 


using  methodologies  incorporating  three-dimensional  vortex-lattice  and  lifting 
line  theories  and  two-dimensional  analog  techniques  which  were  coupled  with 
boundary  layer  prediction  methods  and  empirically-based  ejector  augmentor 
design  and  performance  procedures.  The  objective  of  the  design  was  to 
produce  a  wing  which  would  achieve  high  lift  over  a  range  of  high  angles  of 
attack  in  subsonic  flight  up  to  Mach  nutrtsers  of  0.3  (7). 

Purpose 

As  an  aide  to  assessing  the  data  generated  by  the  testing  of  this  wing, 
it  was  requested  that  a  simple  theoretical  model  of  an  ejector  wing  be 
developed.  This  model,  hopefully,  would  enable  quick  and  easy 
determination  of  the  effects  on  the  lift  of  the  wing  that  result  from 
varying  its  geometric  characteristics.  These  characteristics  include,  for 
example: 

•  the  height  of  the  slot  (i.e.,  the  vertical  distance 
between  the  upper  and  lower  airfoil  sections) 

•  the  horizontal  distance  between  the  leading  edges  of 
the  two  airfoils 

•  the  length  of  the  lower  airfoil  with  respect  to  the 
upper 

The  development  of  this  model  is  the  objective  of  this  thesis.  The 
methodology  of  the  development  and  the  results  obtained  are  discussed  below. 

Geometric  Model 

A  two-dimensional  geometric  model  was  constructed  as  shown  in  Figure  4. 

It  was  based  upon  a  two-element  lifting  surface;  with  a  sink  and  source  to 

represent  the  ejector  inlet  and  exhaust  mass  flows,  respectively  (5,  8,  5). 

It  was  decided  to  simulate  the  upper  and  lower  airfoil  sections  as  parallel 

symmetric  airfoils  and  to  represent  them  by  their  chord  lines.  Each  chord 

line  was  assumed  to  be  of  unit  length,  and  the  leading  edge  of  the  lower 
line  was  selected  to  be  0.55  unit  behind  the  leading  edge  of  the  upper  line. 
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I. 


(Orientation  of  the  two  surface  elements  was  arbitrary  (8).)  This  resulted 
in  an  overall  chord  length  for  the  simulated  ejector  wing  of  1.55  units. 

The  vertical  spacing  between  the  lines  was  assumed  to  be  about  one-tenth  of 
their  collective  horizontal  length  (7).  Thus,  a  vertical  spacing  of  0.15 
unit  was  selected. 

The  horizontal  location  of  the  sink  was  determined  by  placing  it  near 
the  leading  edge  of  the  lower  line;  at  a  distance  of  0.05  unit  in  front  of 
the  lower  line.  The  horizontal  location  of  the  source  was  determined  by 
placing  it  midway  between  the  trailing  edges  of  the  lines.  The  vertical 
locations  of  both  were  established  by  placing  them  midway  between  the  chord 
lines. 

Figure  5  shows  the  complete  geometric  model  of  the  ejector  wing  and 
identifies  all  the  numerical  constants  that  are  involved  in  defining  its 
configuration.  The  configuration  is  drawn  with  respect  to  a  universal  X-Y 
coordinate  system,  but  it  should  be  noted  at  this  time  that  two  other 
coordinate  systems  also  are  present.  They  have  significance  only  in  the 
horizontal  direction  and  they  have  origins  at  the  leading  edges  of  the 
chord  lines.  They  are  designated  the  and  systems  and  will  be 

used  later  to  develop  the  mathematical  model  of  the  wing. 
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LIFTING  SURFACE 


<• 


•  SINK  REPRESENTS  EJECTOR  INLET  FLOW 

•  SOURCE  REPRESENTS  EJECTOR  EXHAUST  FLOW 


Figure  4.  Schematic  of  Ejector  Wing  Mcoel 
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(horizontal  location  of  sink) 

(horizontal  location  of  source) 

(horizontal  location  of  lower  leading  edge) 

(horizontal  location  of  lower  trailing  edge) 

(vertical  location  of  upper  line) 

(vertical  location  of  lower  line) 

(length  of  upper  line) 

(length  of  lower  line) 


Figure  5.  Geometric  Model  of  a  Simulated  Ejector  Wing 


II.  Theory 


General  Equation  Develop'nent 

The  mathematical  model  for  an  ejector  wing  was  developed  from 
consideration  of  the  principles  of  low-speed  aerodynamics.  Central  to  this 
development  was  the  assumption  that  the  viscous  boundary  layer  surrounding 
the  wing  was  thin  and,  therefore,  had  a  negligible  influence  on  the 
inviscid  flow  field  (1).  In  addition  to  ignoring  the  effects  of  viscosity, 
it  was  also  assumed  that  the  flow  about  the  wing  was  steady, 
incompressible,  and  irrotational .  Thus,  for  this  combination  of 
conditions,  the  governing  equation  for  the  flow  field  is  Laplace's 
Equation;  ^ 

where  <^is  the  total  velocity  potential  (1). 

The  model  was  constructed  by  the  addition  of  a  point  sink,  a  point 
source,  and  two  parallel  vortex  sheets  to  a  uniform  stream.  The  vortex 
sheets  are  needed  on  the  chord  lines  in  order  to  produce  a  lift-generating 
pressure  difference  between  the  upper  and  lower  surfaces  of  the  airfoils 
they  represent.  The  vortex  sheets  are  parallel  to  the  uniform  stream  and 
the  sink  and  source  are  located  in  the  area  between  them.  The  procedure 
used  to  develop  the  model  began  with  determining  expressions  for  the 
y-components  of  velocity  (ir)  at  all  points  on  the  sheets.  Since  the 
governing  equation  is  linear,  the  velocity  potential  at  any  point  in  the 
field  resulting  from  this  combination  of  elements  is  simply  the  sum  of  the 
velocity  potentials  for  each  one  (4).  Therefore,  the  velocity  at  any  point 
consists  of  contributions  from  the  individual  components.  Thus,  for  the 
case  involving  no  angle  of  attack,  two  scalar  equations  for  v*  result: 
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"^upper  ®^ink  Source  "’‘^^pper 
sneet 


lower 


vQrti^ity  vprtigity 

distribution  distribution 


(2) 


lower  “  sink  source  ’•'^pper 
sheet  vprtiQity 

distnbuti 


lower 

vprticity 

on  distribution 


(3) 


Considering  first  the  contributions  of  the. sink  and  the  source,  an 


expression  for  the  velocity  potential  at  any  point  (x,  y)  in  the  flow  field 
is: 


^  ®*^sink  ■‘■^source 

(4) 

where: 

4^sink  “  "  -^sink  If'  f“a 

ZTT 

and 

(5) 

4^source  “  -^source  1"  f*b 

(6) 

ZTT 

For  points  on  the  sheets,  expressions  for  rg  and  ri,  in  terms  of  x  and  y  can 
be  written.  They  are: 


Ta  *  ^{x  -  Xj)2  +  y2 

H)  =  (x  -  X2)2  +  y2 


and 


(7) 

(8) 


Substituting  into  the  expression  for<P,  gives: 

<P  •  -Ailnk  in  Jx2  -  2xin  .  xi2  t  y2  * i„ |x2  -  Zxjx  t  xjZ  t  y^’ 

zrr  in 

An  expression  for  '\J'1s  obtained  by  differentiating  this  expression  with 
respect  to  y.  Thus: 


(10) 


IT  ■  "  si  nk 

.y  ■] 

4.  source 

_ y 

aiT 

x2  -  2xix  +  Xi2  +  y2 

fir 

x2  -  2x2X  +  X22  +  y2 

(11) 
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Next  to  be  considered  are  the  velocities  induced  by  the  vortex  sheets 
at  points  on  the  sheets,  themselves.  It  is  known  that  these  velocities  are 
perpendicular  to  the  sheets  (1).  Considering  an  incremental  element  of  the 
upper  vorticity  distribution  (d^i),  the  incremental  velocity  it  induces  at 
a  point  on  the  upper  sheet;  namely  at  (x,  yi)is: 


(12) 


where: 


'i'i(Ji)  =  upper  vorticity  distribution 
=  X  -  5i 

Substituting  for  ri  in  Equation  (12)  and  integrating  over  the  entire  length 


of  the  upper  sheet  leads  to: 

=  fMiL  . 

Wi.)  2rr  X  -  5i  ^^1 


(13) 


Likewise  for  the  lower  sheet,  consider  the  incremental  velocity  induced  by 
6^2  3  point,  (x,  y2).  Integrating  over  the  entire  length  of  the  lower 
sheet  leads  to: 


V, 


J-rTY', 


JITT 


J- 


(14) 


Again,  this  velocity  is  perpendicular  to  the  sheet  and  thus,  is  in  the 
y-direction. 

Additional  complexity  is  found  when  considering  the  velocities  induced 
by  d5i  at  points  on  the  lower  sheet  and  by  djg  at  points  on  the  upper 
sheet.  This  complexity  is  caused  by  the  geometry  of  the  model  but  is 
easily  handled  by  separating  and  focusing  on  the  components  of  the  induced 
velocities  that  are  normal  to  the  sheets,  only. 

Consider  the  velocities  Induced  at  points  on  the  lower  sheet  by  d^i 
(See  Figure  6).  For  this  case,  the  velocity  in  question  is; 
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(15) 


dVy2yd5i  =  4(^l) 
xrrr^ 

where : 

1*3  =  ^(x  -  5i)2  +  (yi  -  y2)2  (16) 


The  angle,  B,  is  involved  in  defining  the  normal  component  of  this  velocity 
and  this  angle  is  defined  as: 


Thus,  the  normal  component  of  velocity  is: 


Integrating  over  the  length  of  the  upper  sheet  gives; 


(18) 


V, 


'  >Tr 


%(l){ tos  [w*  } 


(19) 


>1 


AS. 


Likewise,  the  velocities  induced  at  points  on  the  upper  sheet  by  d^z 
are  (See  Figure  7): 


dVyi/d^o  =  ^2  (^2)d^2 
^  2Trr4 


(20) 


where; 

v“4 « (x  -  X3  -  ^2)^  +  (yi  "  y2)^  (21) 

The  angle  ^  is  involved  in  defining  the  normal  component  and  this  angle  is 
defined  as: 

I  •  tan  -ifyi  -  y2  \ 

-  X3  -  j2  y 

The  normal  component  of  the  incremental  velocity  then  is  seen  to  be: 

dVyi/d J2/n  “  **^2  cos  S 

Integrating  over  the  length  of  the  lower  sheet  gives: 


12 


V* 


^2  (^2)' 

c  * 

'un  -yj'i  -  ^2  'i 

7 

277‘ 

1  ' 

cos  j 

U  -  X3-  ?2l 

i 

-  X3  -?2)^  +  (yi  +  n)^  (22) 

In  order  to  find  the  total  y-components  of  velocity  at  points  on  the 
vortex  sheets,  it  is  necessary  to  sum  the  contributions  from  all  the 
elements  in  the  model  as  per  Equations  (2)  and  (3).  Thus,  for  the  upper 
sheet : 


She^-f- 


1 


iaufct. 


>-rr 


Al _ 

C _ ] 

ATT  , 

c/4 

(23) 


And  for  the  lower  sheet: 


1, 


fouHr 

Skeef 


x^-Xj 


c^4 


(24) 


Equations  23  and  24  are  the  heart  of  the  ejector  wing  model ,  and  the 
method  of  their  solution  is  the  key  to  determining  the  lift  per  unit  span 
of  such  a  wing.  The  method  of  solution  used  involved  satisfying  both  a 
boundary  condition  and  an  auxiliary  condition  and  assuming  algebraic  forms 
for  the  two  unknown  vorticity  distributions:  (Jj)  and  ^2(^2)* 

First,  the  boundary  condition.  In  order  to  make  the  geometric  model 
actually  representative  of  the  aerodynamic  characteristics  of  an  ejector 
wing,  Thin-Airfoil  Theory  states  that  the  camber  lines  (which  for  the  case 
of  symmetric  airfoils  are  coincident  with  the  chord  lines)  must  be 
streamlines  of  the  flow  field.  In  order  to  make  the  cairtber/chord  lines 


14 


streamlines,  it  is  necessary  that  all  velocity  components  normal  to  them  be 
zero  (1).  This  is  the  flow-tangency  condition  and  its  imposition  on 
Equations  (23)  and  (24)  gives: 


The  auxiliary  boundary  condition  is  the  Kutta  condition  which  demands 
that  the  vorticity  distributions  at  the  trailing  edges  of  the  airfoils  be 
zero  (1).  This  condition  was,  satisfied  through  specification  of  the 
coefficients  of  the  functional  expressions  assumed  for  and  X’aiSz). 

It  was  decided  to  use  power  series  representations  for  both  of  these 
functions  and  to  expand  the  series  about  the  trailing  edges  of  the  airfoils. 
Thus,ll(i(^l)  is: 

=  Aq  +  Ai(5j  -  c)  +  A2  di  -  c)2  +  .  .  .  +  An  (Ii  -  c)"  +  .  .  .  (27) 

where  c  is  the  length  of  the  upper  sheet.  Similarly ,^2(^2)  is: 

“  Bq  +  Bi(52  )  +  B2(52  -i)2  +  .  .  .  +  Bn(T2  -/^)'’  +  •  .  .  (28) 

where is  the  length  of  the  lower  sheet. 
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By  assuming  that^l(5l)  and  sf'e  able  to  be  represented  by  power 

series,  it  has  further  been  assumed  that  all  points  on  the  sheets  lie  within  the 
convergence  sets  of  two  series.  The  result  of  these  assumptions  is  that  both 
^l(ri)  and/2(:'2)  are  continuous  over  the  lengths  of  the  sheets.  As  will  be 
shown  later,  these  two  functions  are  indeed  continuous  over  the  lengths  of  the 
sheets  thereby  justifying  the  assumptions  made. 

The  Kutta  condition  demands  that  :^i(c)  =  0  and  ^2.{l)  =  0.  These  conditions 
are  satisfied  by  setting  both  Aq  and  Bq  equal  to  zero.  Thus,  the  series  reduce 
to: 

)(l(?l)  =  Ai(2’i  -  c)  +  A2(?'i  -  c)2  +  .  .  .  +  AnU'i  -  c)^  +  .  .  .  (29) 

and 

'tzi^l)  ~  Bi(^2  B2(^'’2  +  .  .  .  +  Bn(J'2  --^)"  +  •  •  •  (30) 

The  next  step  in  the  solution  of  the  two  major  equations  of  the  model , 
Equations  (25)  and  (26),  involves  substitution  of  Equations  (29)  and  (30)  into 
them  and  the  subsequent  determination  of  Aj  through  Ap  and  Bj  through  Bn  which 
force  the  vortex  sheets  to  become  streamlines  of  the  flow  field.  These 
coefficients  are,  of  course,  infinite  in  number  and,  therefore,  require  an 
infinite  number  of  equations  for  their  complete  determination.  Needless  to  say, 
this  is  an  impossibility  and  truncation  of  the  series  representations  at  some 
finite  power  must  be  done.  It  will  be  shown  later  that  the  accuracy  and 
validity  of  the  lift  estimations  calculated  for  the  model  are  highly  dependent 
on  the  order  of  the  truncations  used  and  that  very  high  order  truncations  are 
required  for  believable  results.  However,  for  the  present  discussion,  it  will 
be  assumed  that  no  truncation  is  used  and  that  an  infinite  number  of 
coefficients  will  be  determined. 

From  inspection  of  Equations  (25)  and  (26),  it  is  seen  that  the  terms 
associated  with  the  source  and  the  sink  take  on  finite  values  for  any  value  of  X 
selected.  (It  should  be  remembered  that  Xj,  X2,  Yj,  ^‘'’‘^-^sink 


16 


are  all  constants  of  the  model  and  are  assigned  finite  values.)  Thus,  in  order 
to  determine  values  for  through  and  Bj  through  Bp,  it  is  necessary  to 
evaluate  the  sink-  and  source- induced  velocities  at  n  different  values  of  X. 
This  is  easily  done  and  the  values  which  result  are  labeled  Nyppgp  and 
Equations  (25)  and  (26)  then  become: 


Since  both  Nipper  ^nd  N]gy^er  can  take  on  an  infinite  number  of  values  (i.e. ,  one 
value  for  each  value  of  X  selected),  it  is  the  generation  of  these  constants 
that  enables  an  infinite  set  of  algabraic  equations  to  be  formed. 

Substitution  of  Equations  (29)  and  (30)  into  Equation  (31)  gives: 


X-S. 


[b.(  ,  (35) 

. 


By  expanding  and  simplifying  the  Integral  expressions  in  this  equation,  the 

parent  equation  for  a  set  of  n  algebraic  equations  is  arrived  at.  Namely,  this 
procedure  leads  to: 
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B,  (%■-() *  S>  l  (Sy/)"iJS^^. . .  *  eda^.j)' 6  JS^  .  V 

1.  ...t*.-;(^l 


'U-yy'^2S-^(y:jS 


Similar  substitution,  expansion,  and  simplification  of  equation  (32) 
produces: 

c,  c.  ^ 

fi,l  dr^)&A  *A>j  a-cTsJi,  *•  •  • "  (i.-cfejs,  * 

x-xrs/-^  y^-yyS/^^  ••* 


Equations  (36)  and  (37)  are  much  more  easily  represented  by  use  of 
matrix  notation.  By  designating  the  integrals  contained  in  these  two 
equations  as  j  and  by  adding  the  subscript,  i,  to  the  Nipper  and  N^ower 
constants,  the  following  matrix  equation  is  arrived  at; 

»1  ,2  •  •  •  Uj  Aj  ^upper^ 

^2 ,1  ^2 ,2  •  •  .  ^2  ,n  ^2  ^pperp 


’upperi 


^n,!  ^n,2  •  •  •  Up^n 


”upper2 

Nuppern 
*^1  owerj 
Nlower2 

^lowern 
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Or,  In  more  compact  notation: 


[U]  [AB]  =  [N]  (39) 

The  solution  of  this  matrix  equation,  namely,  the  determination  of  the 
vector,  AB,  is  accomplished  by  premultiplying  the  N  vector  by  the  inverted 
U  matrix.  Therefore,  the  solution  is: 

[AB]  =  [U]-l  [N]  (40) 


As  stated  earlier,  the  coefficients,  Ai  through  Ap  and  Bj  through  Bp,  which 
comprise  the  AB  vector,  are  the  assumed  coefficients  of  the  two  unknown 
vorticity  functions.  The  values  of  these  coefficients  which  result  from 
this  method  of  solution  are  the  ones  which  cause  the  two  vortex  sheets  to 
be  streamlines  of  the  flow  field.  When  these  coefficients  are  substituted 
into  the  vorticity  functions,  polynominals  of  degree,  n,  result.  The  lift 
per  unit  span  of  the  airfoil  section  is  then  calculated  from  these 
vorticity  functions  according  to  [4]: 


L-  -  pV. 


t  ,) 

^Vi(5i)dJi  +J 


(41) 


Evaluation  of  Integral s. 


The  preceding  section  has  provided  the  general  theory  needed  to 
establish  the  solution  framework  for  the  ejctor  wing,  two-dimensional,  lift 
calculation  problem.  The  solution  of  Equation  (39)  was  shown  to  be  a 
fairly  easy  and  straight-forward  procedure.  However,  the  problem  is 
complicated  by  the  need  to  evaluate  all  of  the  Integrals  shown  In  Equations 
(3b)  and  (37)  and  which  comprise  the  U  matrix  In  Equation  (39).  Basically, 
these  Integrals  are  of  two  different  types  and  are  handled  in  two  different 
ways.  The  first  type  is  associated  with  the  Ai  through  A^  coefficients  In 
Equation  (36)  and  the  B^  through  Bp  coefficients  In  Equation  (37).  It 


Is  of  the  general  form: 
rl- 


f(M 

a 


d5 


•'I 

and  possesses  a  singularity  at  the  point,  I  =  a.  In  attempting  to  evaluate 
the  integrals  of  this  type  numerically,  the  presence  of  the  singularity 
created  somewhat  of  a  problem. 

The  second  type  of  integral,  although  more  complex  in  appearance,  does 
not  contain  a  singularity  and  thus,  is  much  more  amenable  to  numerical 
methods.  This  type  is  associated  with  the  coefficients  through 
Equation  (36)  and  the  coefficients  through  Ap  in  Equation  (37).  It  is 
of  the  general  form: 


fil)  9(5)  dS. 


Regarding  the  first  type,  it  was  decided  to  evaluate  those  integrals 
analytically.  The  most  complicated  of  these  is,  of  course,  the  one  of  nth 
order  and  its  evaluation  will  be  presented  as  representative  of  all  others. 
The  nth  order  integral  is  of  the  form: 

,a- 


(5-  c)n 
a  -5 


Expansion  of  the  numerator,  for  all  values  of  n,  produces  polynomials  that 
have  coefficients  which  are  determined  according  to  Paschal's  Triangle. 
Ignoring  these  coefficients  for  the  time  being,  this  expansion  shows  that 
the  general  form  of  the  integral  is  actually  made  up  of  n+1  individual 
Integrals  as  follows: 
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From  this,  it  is  seen  that  the  original  problem  of  evaluating  one  integral 
containing  (5-  c)^  reduces  to  evaluation  of  n+1  integrals;  each 
containing  one  term.  Fortunately,  an  analytical  expression  is  available  to 
aid  in  evaluating  the  single  integrals.  Their  general  expression  is  of  the 
following  form: 

where: 

b  =  -1 
and  m  =  1. 


The  analytical  expression  used  to  evaluate  these  integrals  is  [2]: 


[a  +  b^r 


nU-a)S(a  b5)n“m-s-H 
(n  -  s)!s/(n  -  m  -  s  +  1) 


except  when  n-m-s  +  l  =  0.  In  that  case,  the  corresponding  term  in  the 
square  brackets  is: 


(n  -  m  +  1) !  (m  -  1) ( 


a  + 


With  the  aid  of  these  formulas,  an  integral  containing  any  power  of 
can  be  evaluated.  The  expressions  Tor  the  integrals  containing  terms  up  to 
are  shown  in  Table  I. 

The  second  type  of  Integral  was  evaluated  by  way  of  a  numerical 
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Expansions  of  Integrals  of  the  Form 


integration  method  commonly  known  as  the  trapezoidal  rule.  This  method  is 
used  for  functions  which  are  continuous  over  an  interval  and  involves 
evaluation  of  the  function  at  points  within  the  interval.  Given  the  values 
of  the  function  at  the  points  selected;  namely  fo.  ^1#  ^2.  •  •  .  fn*  the 
trapezoidal  rule  gives  the  value  of  the  integral  according  to  the  following 
expression  [3]: 


f(x)dx- 


^0  '''  ^^2  *  •  *  '*'  ^^n+1  ^ 


nj  2 


(43) 


where  Ax  is  the  spacing  between  evaluation  points. 

The  two  procedures  described  above  provide  all  the  information  needed 
to  evaluate  integrals  encountered  in  attempting  to  solve  this  problem.  All 
that  remains  is  to  determine  how  many  integrals  there  will  be.  The  answer 
to  this  question  is  directly  dependent  on  the  order  of  truncation  of  the 
power  series  used  to  represent  the  vorticity  functions  and  the  number  of 
integrals  increases  rapidly  with  increasing  truncation  order.  It  turns  out 
that  for  truncation  at  order  n,  (2n)2  integrals  are  required.  In  order  to 
generate  these  integrals,  n  separate  control  points  need  to  be  selected. 

Control  points  can  be  any  points  on  the  vortex  sheets  except  their 
beginning  and  end  points.  The  use  of  either  of  these  two  points  gives  rise 
to  integrals  which  are  undefined  at  either  the  upper  or  lower  limits  of 
integration.  These  integrals  do  not  have  finite  values  and,  therefore,  can 
not  be  used  in  the  purely  numerical  computation  procedure  used  to  find  the 
unknown  vorticity  function  coefficients. 
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In  order  to  prove  the  validity  of  the  theory  presented  in  the  preceding 
section,  it  was  decided  to  perform  an  actual  lift  calculation  using  five 
control  points.  This  required  creation  of  the  FORTRAN  computer  program 
contained  in  Appendix  A  and  the  truncation  of  the  assumed  power  series 
representations  of  the  vorticity  functions  at  the  fifth  power.  The  control 
points  were  arbitrarily  chosen  to  be  .2,  .5,  .8,  1.1,  and  1.4  and  were 
selected  so  as  to  fairly  evenly  cover  the  1.55  unit  length  of  the  total 
airfoil  section. 

Before  beginning  discussion  of  this  procedure,  it  is  necessary  to  first 
specify  the  values  of  the  numerical  constants  which  are  required.  The 
geometric  constants  are  the  same  as  those  discussed  earlier  and  are  found 
in  Figure  5.  The  density  and  velocity  of  the  free  stream  were  determined 
according  to  Mach  0.3  flow  at  sea  level.  Thus, 

P  =  .0023769  slug/ft3 
and  Va,=  334.8  ft/sec 

The  values  of  the  sink  and  source  strengths  were  determined  by  considering 
the  two-dimensional  area  rate  of  flow  incident  on  a  line  located  between 
the  vortex  sheets  and  perpendicular  to  the  direction  of  the  free  stream. 
This  area  rate  of  flow  is; 

A  »  (yi  -  ya)  V-  (44) 

and  carries  the  units  of  ft^/sec.  The  sink  strength  was  arbitrarily  set 
at  0.7A.  For  the  source,  it  was  decided  that  its  strength  should  be  six 
times  greater  than  that  of  the  sink.  The  value  of  the  source  strength  also 
Included  the  mass  which  was  "lost"  into  the  sink  and  was  determined 


according  to  the  following  formula: 
-^source  *  J- sink  +  6 sink 


(45) 
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The  reason  for  adding  the  sink  strength  to  that  of  the  source  was  caused  by 
the  need  to  "regain"  the  flow  which  was  absorbed  by  the  sink.  In  reality, 
this  flow  is  associated  with  the  free  stream  and  represents  the  air  which 
passes  through  the  wing  slots.  The  reason  for  using  the  factor  of  six  in 
the  other  source  term  was  to  account  for  the  substantial  mass  addition  to 
the  slot  flow  caused  by  the  ejector.  The  results  of  these  calculations 
are: 

-i-sink  “  35.1540  ft^/sec 
A.  source  “  246.078  ft^/sec 

In  order  to  generate  the  system  of  equations  needed  to  solve  the 
problem  as  posed.  Equations  (36)  and  (37)  with  fifth  order  truncations  are 
used.  These  equations  become: 


x-i;  x-r, 

*A,iys.-c)'eJs,  ♦ 

T?  _  f'l  fr  .P\^ 


(47) 

In  order  to  determine  unique  values  for  through  A5  and  through 
®5,  each  of  these  equations  was  evaluated  for  the  five  different  values  of 
X  corresponding  to  the  control  points  selected.  For  example,  for  the  first 
control  point,  X  ■  0.2,  the  two  equations  became: 
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+  5,  j\s.-^)aJi^  *  B,jf(s,-j)‘tJs^->- 

B^I\s,-/Y6,JS,  ■<■  B^jf  (!,■/)''& + 

A  ^^(K-<)dJS,*/)^l‘(S.-4sJS,*A3i‘'(^r‘yeJK  * Av  f,(S.-^)''^K- 

Similar  equations  were  generated  for  the  four  other  control  points  by 
replacing  the  constant  0.2  by.  successively,  0.5,  0.8.  1.1,  and  1.4.  From 
these  10  equations,  100  different  integrals  result.  For  example,  for  the 
^l  coefficient,  the  integrals  are: 


Ul,l 


.x-r, 


“  r 


r 

u  .  f*  (^i'^) 

^j  (-5--  r.y  (3,-3*)* 


(50) 

(51) 

(52) 

(53) 

(54) 

(55) 
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The  values  of  these  integrals  comprise  the  first  column  of  the  U  matrix. 
Similarly,  the  integrals  associated  with  the  B5  coefficient  make  up  the 
tenth  column  of  the  matrix.  They  are: 


(53) 
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'li 


* 


Us.io  * 


COs(i».^''  aV-X, -Ja) 


JJ, 


U6 


“7.10  ’X 


■  1 77:'^^ 
.10  -\ 


Ulo 


(64) 


(65: 


(66] 

(67: 

(68: 

(69] 


The  other  80  integrals  are  defined  in  similar  fashion  for  the  other 
coefficients,  but  they  will  not  be  listed  here.  It  is  sufficient  to  say 
that  the  U  matrix  which  results  from  these  integrals  is  10  by  10  in  size 
and  can  easily  be  inverted  by  existing  routines  (3). 

The  values  of  Nypper  and  N^Qy^gr  were  calculated  for  the  five  different 
control  points  from  Equations  (33)  and  (34).  These  constants  comprise  the 
N  vector  in  Equation  (39).  From  Equation  (40),  values  for  the  unknown 
coefficients  were  determined  to  be; 

Ai  -  -130782.2 
A2  »  -4213634.4 
A3  «  -20230426.5 
A4  -  -32400910.7 
As  »  -16521998.5 
28 


Bj  ■  -45936.3 
B2  =  545461.9 
B3  =  4372389.1 
B4  =  8572411.4 
B5  =  4925628.4 

These  coefficients  give  rise  to  vorticity  functions  that  cause  the 
y-components  of  velocity  at  the  control  points  to  be  zero.  The  lift  per 
unit  span  is  calculated  from  these  functions  according  to  Equation  (41)  and 
is: 

L'  =  -2230.5  Ibf/ft. 

The  negative  sign  associated  with  this  value,  needless  to  say,  looks  bad 
and  casts  doubt  on  the  validity  of  the  procedure  used.  However,  inspection 
of  Figures  8  through  11  shows  that  the  boundary  and  auxiliary  conditions 
have  been  satisfied.  Namely,  the  y-components  of  velocity  are  zero  at  the 
control  points  and  both"^!  (^i)  and  ^2 (-2)  are  zero  at  the  trailing  edges  of 
the  vortex  sheets.  Thus,  in  spite  of  the  seemingly  wrong  answer  for  L',  it 
appears  that  the  model  and  theory  are  correct  and  that  the  computer  program 
performs  as  expected.  As  a  step  toward  investigating  potential  reasons  why 
negative  lift  was  obtained,  the  calculations  were  repeated  for  different 
selections  of  control  points.  The  results  obtained  are  summarized  in  Table 
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Figure  9,  Lower  Vorticity  Distribution  for  Control  Points:  0,2,  0.5,  0.8 


TABLE  II 


Lift  Calculations  Resulting  from  Four  Different  Control 
Point  Selections  for  Fifth  Order  Truncations 


Control  Points 

,2,  .5,  .8,  1.1,  1.4 
,1,  .45,  .80,  1.15,  1.50 
.05,  .40,  .75,  1.10,  1.45 
.1,  *T,  .9 


Lift  per  Unit  Span 
-2230.5  Ibf/ft 
3690.1  Ibf/ft 
-289.3  Ibf/ft 
343.8  Ibf/ft 


From  Table  II,  it  is  apparent  that  the  lift  is  directly  related  to  the 
control  points  used.  Of  course,  there  are  an  infinite  number  of  possible 
control  point  selections  and  the  criterion  by  which  the  “correct"  ones 
could  be  determined  is  unknown.  Further,  since  the  choice  of  control 
points  was  arbitrary,  their  selection  should  not  have  an  impact  on  the 
final  answer.  However,  the  consistency  with  which  the  boundary  and 
auxiliary  conditions  were  satisfied  does  appear  to  indicate  that  the 
solution  method,  itself,  is  valid  (See  Figures  12  through  23). 
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Figure  18. 
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F-fgure  19.  Y-Components  of  Velocity  over  Lower  Sheet  for  Control  Points: 
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IV.  Conclusions 

A  conclusion  that  can  be  drawn  from  consideration  of  the  numerical 
examples  is  that  more  than  five  control  points  are  needed  for  determination 
of  the  lift.  It  has  been  shown  that  the  flow  tangency  boundary  condition 
is  satisfied  at  the  control  points,  but  it  is  not  satisfied  at  other  points 
on  the  sheets.  It  is  believed  the  validity  of  the  lift  calculated  from  the 
model  is  directly  dependent  on  how  well  the  sheets  are  made  to  coincide 
with  streamlines  of  the  flow  field.  From  inspection  of  the  plots  of  the 
y-components  of  velocity  over  the  sheets,  it  is  apparent  that  for  the 
five-control -point  case,  the  sheets  do  not  closely  resemble  streamlines. 

By  using  more  control  points,  the  flow  tangency  requirement  would  be 
satisfied  at  more  points  on  the  sheets  and,  as  a  result,  they  would  more 
closely  resemble  the  streamlines.  Use  of  more  control  points  will 
entail  much  more  work,  however.  This  work  will  not  be  theoretical  in 
nature,  but  rather,  will  involve  the  fairly  tedious  procedure  of 
determining  values  for  the  definite  integrals, 

r  js 

\  a-j 

where  n  takes  on  values  from  zero  to  a  very  high  integer. 
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V.  Reconmendations 


t. 


Inspection  of  Table  I  shows  that  the  series  which  result  from  the 
integrals  of  concern  are  not  built  one  upon  the  other  —  i.e. ,  by  simply 
adding  a  new  term  to  a  previous  string  of  terms.  It  is  seen  that  the 
coefficients  for  any  term  containing  (a  +  b^)  to  any  power  change  as  n 
changes.  This  lack  of  simplicity  in  the  generation  of  these  series 
complicates  the  manner  in  which  they  could  be  developed  by  a  computer.  If 
these  series  and  the  associated  integrals  are  worked  out  by  hand,  the  large 
number  of  terms  involved  can  soon  become  overwhelming.  For  example,  an 
integral  containing  only  a  fifth  order  polynomial  gives  rise  to  six 
individual  integrals  containing  the  following  terms:  x^,  x^,  x^,  x^,  x^, 
and  x®.  These  integrals  contain,  respectively,  six,  five,  four,  three, 
two,  and  one  term  and,  the  original  integral  requires  the  sumnation  of 
these  terms  —  i.e.,  21  —  for  its  evaluation.  Clearly,  it  can  be  seen 
that  to  do  this  procedure  by  hand  would  be  very  laborious. 

The  most  promising  method  for  solution  of  this  problem  appears  to 
involve  creation  of  a  computer  routine  based  on  Equation  (42)  which  would 
generate  the  series  representations  of  the  integrals  for  any  value  of  n. 

If  such  a  routine  were  available,  many  more  than  five  control  points  could 
be  used  and  thus,  the  vortex  sheets  could  be  brought  more  into  line  with 
the  flow  field  streamlines.  Then,  the  calculated  values  of  L*  probably 
would  not  display  the  same  erratic  behavior  as  they  have  done  for  the  cases 
where  only  five  control  points  were  used. 
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APPENDIX  A 


Listed  on  the  following  11  pages  Is  the  computer  program  that  was  used 
to  make  the  lift  calculation.  The  Inputs  to  the  program  are  contained  In 
lines  17  through  39  and  consist  of  the  required  geometric  values  and 
aerodynamic  constants.  Evaluation  of  integrals  is  done  in  lines  44  through 
252,  and  the  sink-  and  source-related  constants  are  generated  in  lines  260 
through  273.  Following  these  calculations,  the  U  matrix  is  inverted  by 
means  of  the  subroutine  labeled  INVDET. 

Output  from  the  program  begins  with  line  315  which  prints  out  values  of 
the  upper  vorticity  distribution.  Similarly,  line  344  prints  out  values  of 
the  lower  vorticity  distribution.  The  overall  lift  per  unit  span  is  output 
via  line  364. 

The  program  also  contains  a  large  section  which  verifies  that  the  flow 
tangency  condition  is.  Indeed,  satisfied.  Velocities  and  pressures  at 
points  along  the  upper  sheet  are  output  via  line  459  and  similar  values  for 
the  lower  sheet  are  output  via  line  538. 

The  plots  included  in  this  report  were  not  generated  by  this  program, 
but  were  made  via  a  separate  plotting  program.  Potential  users  of  this 
program  are  advised  to  procure  similar  independent  plotting  routines  If 
graphical  output  Is  desired. 
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C 

C 

CALL  INVOFTlU. ln.0TNRM.0ETM) 

WR1TCI8.210) 

210  FORMAT!*  THE  INVERTEO  0  MATRIX  IS!*) 

00  220  I ■  1.10 
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UNCLASSIFIED 


28S 

09ITE(0.*>  (UII.J),  J-I.IO) 

7a« 

220 

CONTINUE 

207 

c 

200 

c 

209 

00  230  I>l  .  io 

290 

AB(j|aU(I.t>«N(l)«U(l.2l«N(2)^U(l,3)*Nt3)^ 

29| 

C  U(t,9l*N(9)>U(I.5}«N(5>^U(Ii0)*N(0) 

292 

C  ♦UI 1 ,7J*NI7)*U( I ,8)«N(8) 

293 

C  ( I  I  9 ) •N ( 9 1 «UI I • 1 0 ) ^N ( 1 0 ) 

299 

230 

continue 

295 

WRITE(0t29OI 

29t> 

290 

FORMAT!*  THE  COeFFS  .  A  1  •  A2  .  A3  ■  A9  ,  AS • B  1  .  02 i B3 • 89 . b5  ARe:*/» 

297 

00  250  I » 1 . 1 0 

290 

*R1TF<  *1  •  (  1  X  ,F2a,  a  >  *  1  ABin 

299 

250 

continue 

3oq 

c 

301 

c 

NO*  FOR  The  lift  calculation 

302 

c 

303 

NLlFTalOO.O 

309 

STEPL«100 

305 

SUMbO.O 

300 

TeHP»0.0 

307 

NUMRRbO 

300 

• 

X-0.0 

309 

INCL1«(C-0.0)/N|.IFT 

310 

00  200  I»1 ,STEPL+1 

31  i 

GAHMAl»Aa<  U*<*-C)*A8(2J«(CX-C>**2»+AB<3»*( (X-C>**3) 

312 

C  ♦AB{9J*(  (X-C)**9)*A6(5>*( 

313 

c 

319 

c 

315 

*RITE(3.2701  X.SAMMAl 

310 

270 

FORMAT! IX,2Fl0t8) 

317 

c 

310 

c 

319 

NUMaRsNUMbR^I 

320 

IF  I!NUMBR  .EO*  11  (OR*  INUHBR  tEO*  STEPL^I))  THEN 

321 

SUMsSUM^GAMMA 1 

322 

else 

323 

tehp«gammai*2.o 

329 

SUM«SUM  +  TeH«» 

325 

ENO  IF 

320 

XaX^lNcLI 

327 

200 

CONTINUE 

.320 

Ll>SUM«! INCLt/2. I 

329 

*RITE!0«28O) 

330 

200 

FORMAT!*  THE  VORTICITY  OF  THE  TOP  PLATE  IS!*/) 

331 

*(R|TCI0t*l/lX,F2O.IO)»)  Ll 

332 

SUMaO.O 

333 

TFMp-O.O 

339 

NUNpRaO 

335 

c 

330 

XaO.O 

337 

INCL2-Ix9-X3>/NLIFT 

330 

00  290  I**ltSTEPL*l 

339 

6ANMA2aAR!0)«|X.BLl*AB(7)«(X-BLI**2«AB(8)*(X*BL)**3 

390 

C  ♦AB(9)«Ix«RLI**9«AB(10)*(X-BLI**5 

39i 

c 
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uricLASsirten 


3*t5 

3**6 

3M7 

S'tet 

3*«9 

350 

351 

352 

353 
359 

355 

356 

357 

358 

359 
380 
3<^t 
3A2 
3«3 
389 
385 
388 

387 

388 

389 

370 
37  I 

372 

373 
379 
375 
378 

377 

378 

379 

380 

381 

382 

383 
389 
385 
388 

387 

388 

389 

390 

391 

392 

393 
399 
395 
398 

397 

398 


C 

C 

8RITE<  1  I  *3C0»  *,SaMMA2 
300  FORMAT! IX. 2Fl8*a) 

C  . 

C 

NUMRRSNUM8R4. 1 

IF  ((NUHBR  .Ed*  1)  tOR*  (NUH8R  .£0.  STEPL^-t))  ThEn 

SUMaSUM«G AMMA2 

ELSc 

TEHPaGARRA2«2>Q 

SUH«SUm*TeMP 

end  if 

X»X* iNcL2 
290  CONTINUE 

L2«SUM* ( INCL2/2.  ) 

WRITC(8.310) 

310  FORMAT!*  THf  VORTiCITY  OF  THE  BOTTOM  PUATE  ISt*/) 

RRITE ! 8, •! /fx ,F20. 10) • I  L2 
CIFT«i!RH0*VeL)*(1.1*L2) 

»RITE!8.320) 

320  FORMAT!*  THF  VALUE  OF  THE  LIFT  ISt*/) 

WRITE  !  •  .  *  I  1  X  .F  I  8..a  )  *  )  LIFT 
C 

c  This  is  the  check  portion  of  the  program 
c 

C1»!-I.)*AB!5) 

C2«A8(9)-5»«C«Aa(5l 

C3*A8(3»-9»*C*A8(9>+1C«*C**2*ABI5) 

C9*A8(2)*3»*C*AR(3)t8»*C*»2*AB(9)«10»*C**3*ABI5) 

C5"AB!i)-2»*C*A8!2>*3«*C**2*ARI3)« 

C  5.*C»*3*AB! 9)*5«*C**9«AB(5) 
C8"C*Ae(l)-C*«2*A8<2)‘*C«*3*AB!3)-C»«9«AB(<|)4' 

C  r**5*An!5) 

C7« ! -1.)*AB(10)_ 

C8«ABI9)-5««BL*aB!10I 

C9>AB(8)-9t«BL«AB(9)^10,*BL*«2*AH(l0l 

C10«A8I7)-3.*9L*AB!8)t8»*BL**2*A8C9)-10»*BL»*3«aB!10) 

CllaAB!8)-2<«BL*A8(7lT3«*BL**2*AB(R)> 

C  9  * ^BLaaB* AB ( 9 ) ^5  •  •8L**9» AB  !  1  0 > 

Cl2«BL*AB!8)-BL«a2*A8!7>aBL«*3aAB!B)-8L*»9«AB!9)a 
C  RL««5«A8(ia) 

C 

C  FOR  the  top  plate  we  have; 

c 

WRITE(8«50C) 

SOO  FORMAT!*  V.U  ANO  PRESSURE  ALONG  THE  TOP  PLATE*/) 

XaO.05 

DO  321  l>tf2G 

Pill  )»SOURCF*  (Yi/(X**2-2,«X2*XaX2»9  2'*Yl*»2>) 

C 

P9I I )•( SOURCE/ I2**P1 ) )*l CX-X2)/IX»92-2*9X2*XaX2**2aY 1**2) ) 

C 

P2li)*»(-l«)«(SlNK*lYi/IX*«2-2»»Xl*XaXi99  2aYl**2)>) 

C 

Pl0(lI*<(-l*«SlNK)/(2«9PI))*((X-X|)/(X**2>2«*Xl*XaXiaa24-yi**2)) 


57 


UNCLASSIFIED 


wx^un-'.ai’'  ICU 


399 

900 

901 
90? 
903 
909 
909 
90A 
907 
909 

909 

910 
91  i 

913 
9  m 

915 

916 

917 

919 
9  1  9 

920 

921 

922 

923 
929 

925 

926 

927 
926 

929 

930 

931 

932 

933 
939 

935 

936 

937 
936 
939 
990 
99  1 

.992 

993 

999 

995 

996 

997 
996 
999 

950 

951 
962 
953 
959 
955 


63<I»«C1*1C»*5/5«>(X*C*»9)/9«*(X**?*C»«3J/3»*(X**3*C**2)/2»* 

C  X«*9*C^X« •S* AL06 ( ABS  I  ( X-C  1 /X  1  )  I - 

C  C2»(C*»9/9«4<X*C»«3)/3*-»<X«*2«C««2>/2»* 

C  X**3*C+X*«9*AL0G(A8StlX-C)/X>)»- 

C  C3«(C**3/3.*{X«C**2)/2.*X**2«C* 

C  X**3*AL0G{ABS<(X-C)/X)n-C9*(C«*2/2#> 

C  X*C4.X**2*AL0GCABSI(X-C>/X>>)-C5*1X* 

C  ALOGl aBS ( < X-C 1 /X) » ♦€! ♦C6»AL0G 1 ABSC ( X-C) /X ) ) 

C 

SUHaOtO 

TfMpaO.O 

Xl2aO.Q 

c 

no  322  K*  1  •  1 0  1 

T0PB(AB(6)a(Xl2-8L»>AB(7)a(XI2-BL)**2* 

C  AB*8|*(X12-BL)**3^AB19)*1XI2-8L)*»9* 

C  AB<IQ)«(XI2-BL»**5)«(cOS(ATAN(H/(X-X3-XI2M)) 

BOTaSORTl  (X*><Xl2'*^X3)  )*«2*H**2> 

01 VaTOP/BOT 
C 

IFllK  .fQ*  1)  (OR*  (K  tCQ*  lOin  THEN 

SllHaSUM«0  I  V 

ELSE 

TEMP«0IV*2« 

SUMaSUM^TCMp 
ENO  IE 

c 

XI2«XI2-»«0l 
322  continue 

C 

P9« i )»SUMa( .01/2*  > 

C 

SUMaO.O 

TcHPaO.O 

X|2a0*0 

c 

00  329  1  I  1 C  1 

TOPal  AB  (6)  •(  X  I2.BL  )«AB  (  7)  *1  X  1  2-BL  1  a*24- 
C  AB181*(XI2-BL»**3*AB(91*(XI2-BL)**9* 

C  AB(10)*(XI2-3U)*«5J«(SIN(ATAN1H/(X-X3-XI2)))1 

BOTaSORTl  (X>(XI2'^X3)  )«*2*H**2> 

OlVaTOP/BOT 

C 

IFIIK  .fQ*  1)  tOR*  (K  «EQ«  1011)  THEN 
SUHaSUM<»OlV 

else 

TEMP»01V«2. 

5UM«SUM«TEMP 
E*<0  IF 

c  . 

XI2>Xt2-»«01 
329  continue 

C 

P|icn»<l./(2.*PI)  l«SUNa(  .01^2,  ) 

c 

VTOPl  I  )>pi  (  T  )4.P2(  I  )9P3(  I  I^PRl  I  } 

UTOPl I )aVFL4P9( i I♦P10( I mP| I (  I  I 
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unclassified 


47 

48 


43 

46 


48 


^0 

n 

72 

n 

78 

75 

76 

77 

78 

79 


96 

97 

98 

99 


^7 


I  I 
12 


C 

PTOPI1)«(.5«RhO«VEL**?»*(1»-((UTOP<I»**2*VTOP(I)**2>/ 

C  « VEL*«2 ) J » ♦PPeS 

write*  6. 3  73  »  X , VTOP* I ) tUTOPi I ) iPTOPC 1  ) 

XaX^.OS 
32!  continue 
C 

C  FOR  THE  BOTTOM  PLATE  «E  HAVE! 

C 

WRITC(6.6G0) 

600  FORMAT!*  ViU  ANo  PRESSURE  ALONG  THE  BOTTOM  PLATE*/* 

X>0.6 

Do  330  l»t.?0 

PSII)aSnURCr«(Y2/(X««2-2.aX2*X>X2«*2^Y2»*2)) 

C 

Pl2Cl!"<S0URCE/<2.*Pl»)*MX-X2>/<X**2-2»*X2aXaX2«*2aY2**2)> 

C 

P6«IJ«(-1.»*(SINK*(Y2/<X**2-2«*X1*X*X1**2*Y2**2»>) 

C 

Pl3(I)«((-l,«SlNKl/(2**PI))*((X*Xl)/(X**2-2»*Xl»XaXl**2aY2**2)> 
C  . 

P8<II«C7*(BL**5/5.'*-<(X-X3)»BL«*8)/8»*'((X-X3)**2*SL**3)/3* 

C  ♦  (  (  X.X3  I  ••3*BLa«2  I /2.*' 

C  (X-X3»»«H*0L*'«X-X3J**5*AL0G(ABSr(X-X3-BL)/*X-X3))))- 

C  Ca*<8L**RA.>((X-X3)*BL**3)/3.  +  («X-X3)**2*BL**2)/2.a 

C  IX-X3>**3»BL*-(x-X3>«*8*AL06(ABS(  (  X  -  X  3-BL  )  /  (  X  “X  3  )  )  )  )- 

C  C9*|BL**3/3,a({X-X3J*BL»*2)/2.'f(X-X3)**2*BL  + 

C  <X-X319a3*AL0G*AflS<  » X-X3-eL ) / ( X-X3 > ) > ) -C ! C* f S L • • 2/ 2 . ♦ 

C  IX-X3>*eL+(X-X3»**2«AL0GU8S(*X-X3-BL)/(X-X3))>» 

C  -Cl !•( (X-*31*AL0G(A8S( (X-X3-BL)/(X-X3)>) 

C  aRLlaCU^ALOGl  ABSl  (X-X3-BL>/(X-X3))) 

C 

SUMaOaO 
TEMPaO.O 
X  I  I a0*0 

c 

00  331  X>1 • 101 

TOPalAB! n*(XIl.C)^AB(2)*(Xll.C)«*2a 
C  AB<31*(XI1-CI**3*A8<H>*<X11-C)**8^ 

C  AB(5I*(XI 1-CI**5»*(C0SIATAN(h/<X-XI  1  )  1  )  ) 

B0TaS0RT!<x-Xlll**2*H*a2) 

OlVaTOP/BOT 

C 

IF((K  tEO*  1)  (OR*  (K  tEQ*  101>)  THEN 
SUMaSUMaOI V 

else 

TEMP**0IV*2« 

SUMaSUMaTEMP 
END  IF 
C 

XI 1>X1 14*01 
331  continue 
c 

P7( I laSUMal .01/2*  1 
C 

SUMaO.O 

TERPaO.O 

59 


unclassified 


•  •  •  W 


SI3 

*11-0.0 

SIH 

C 

SIS 

00  332  tc-1  1  lOl 

SI& 

T0P-( AB( l)»t*ll-C)^ABC2)««XIl-C)**2t 

S|7 

c 

AB«3>*<*ll-Cl**3^AB(9M(xll-C>*»9* 

5U 

c 

ABI5)«(*11-C)»*5>*(S1N(ATAN(H/(X«XI1)))> 

519 

B0T»SQRT<(X-xril**2'»H*-21 

520 

OIV-TOP/BOT 

521 

C 

522 

1FI(K  .PQ.  il  .OR.  IK  .eO.  101)1  then 

523 

SUMaSUM-.01  V 

529 

ELSE 

525 

TEMP»0IV*2. 

52« 

Sum-Sun.temp 

527 

END  IF 

526 

c 

529 

xii-xiu.Gi 

530 

332 

CONTINUc 

531 

c 

532 

Pl9in*(l./<2«*P1)  1  •SUM*  (  .01/2*  ) 

533 

c 

539 

VaOTlIl-PSfp^PAnl^PTlD^PSll) 

535 

UB0T<l)-VFL*P12(l)>P13(l)-*-P19(t) 

S3* 

PR0T*l)«(.S*RH0*VEL*-2)*<l.-(C>J80T(I)*-2*VB0Tcn**2)/ 

537 

C 

( VCL**2l  1  1  ♦PRES 

538 

WRITE!*. 323)  X.VBOTd)  ,UBOT(  I  )  ,PbOT(  1  ) 

539 

X-X^.05 

590 

330 

continue 

591 

c 

592 

323 

EORMATI  IX,9F1A.5) 

593 

END 

59** 

subroutine  iNVOET(C.N,OTNRM,OETM 1 

595 

dimension  CCN.N),  J(10Q) 

59* 

P0"1. 

597 

DO  129  L*I.N 

596 

OO-Q. 

599 

00  123  K- 1 .N 

550 

123 

00-00  +  C(L,IC)*C(L.K) 

551 

OO-SORT(OO) 

552 

129 

PO-PO-OO 

553 

DeTm-1, 

559 

on  1 25  L*l (N 

555 

125 

J(L^20)-(. 

55* 

00  199  L-I.N 

557 

cc-n. 

556 

M-L_ 

559 

00  135  K-L,N 

5*0 

IP  M  ARS(CC)*ABS  (C  (I..K  )  n  .GC.O*  )  60  TO  135 

5*1 

12* 

M-K 

5*2 

CC*C(L,K) 

5*3 

135 

continue 

5*9 

127 

IP  (L.EO.M)  GO  TO  138 

5*5 

126 

K-JI M^20 ) 

SA* 

J(M^2Q)-J(L^20I 

5*7 

J(L420)-K 

5*6 

Oo  137  K-1,k 

5*9 

SaC  I  K  .L) 

60 


uNct«ssf rten 


U'«LI.  A  3  3  I  t  ltd 


570 

571 

137 

C(K,M)aS 

572 

138 

C»L.L)»1» 

573 

oetm-oetm.cc 

57‘« 

00  1 39  M«1  ,N 

575 

139 

CCL,M)-C<LtMl/CC 

57& 

00  192  M*! iN 

577 

If  60  TO 

192 

57e 

129 

CC-C«M.L  > 

579 

IF  (CC*eQ*0«>  60 

TO  192 

580 

130 

C(H.L)»0« 

5At 

DO  191  K«1,N 

582 

191 

C1M.K)«C<M,IC1-CC*C«L*K) 

5«3 

192 

continue 

589 

1  99 

continue 

585 

00  193  L«1 ,N 

586 

IF  ( J(L«20I *EQ*L> 

GO  TO 

587 

131 

MaL 

588 

132 

M«M*  1 

589 

IF  ( J(M-»20l  •eQ*L) 

GO  TO 

590 

136 

IF  (N*6T*M)  go  to 

132 

591 

133 

J(M*20) aj(L^20) 

692 

00  163  iN 

593 

CC»C(U.K) 

599 

C(L,KI«C(N,K) 

595 

163 

C (H ,K 1 aCC 

596 

J(L<»2Q|aL 

597 

193 

continue 

598 

OFTMaABS(OETM) 

599 

0TNRM«0ETM/P0 

600 

return 

601 

END 
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